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ABSTRACT 


The temperature distribution T( r , t) in homogeneous opaque 
spheres with constant radius r = R and constant initial temperature 
T 0 = T(r, 0) is represented as a closed-form solution of Fourier’s 
equation when heat is radiated from the surface r = R into a vacuum. 
Substitution of the boundary conditions at r = 0 and r = R into this 
solution yields a system of two nonlinear integro- differential equations 
for T(R, t) and T(0, t ). These equations are solved numerically. 
Results for T(r, t) are presented for spheres made of quartz with 
R = 1 cm, R = 15.3 cm, and R - One application of the solution is 
to determine the temperature history of meteorite ejecta traveling in 
space. 
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THEORETICAL SOLUTION OF THE NONLINEAR PROBLEM 
OF TRANSIENT COOLING OF AN 
OPAQUE SPHERE IN SPACE 

By 

Ernst W. Adams* and H. Spreuer 
Deutsche Ver sucks anstalt fur Luft- und Raumfahrt 
Freiburg i. Br., West Germany 


INTRODUCTION 

The transient cooling process of spheres with constant radius r and constant material proper- 
ties is considered under conditions of radiative heat transfer from the surface to the environment 
and conduction of heat in the material, which is assumed to be homogeneous and opaque to thermal 
radiation. The particular example which induced this study concerns the temperature history of 
debris ejected from the moon by a meteorite impact. The spherically symmetric closed- form 
solution of Fourier’s heat conduction equation which is employed here is well known for the case 
of linear boundary conditions. The mathematical problem is due to the nonlinear boundary condi- 
tion at r = R, following from the Stefan- Boltzmann radiation law. By making minor changes in the 
derived solution, other forms of this boundary condition may be employed, such as those due to 
additional convective heat transfer and/or mass transfer. In a forthcoming publication, a general- 
ization of the solution presented here will be employed to solve the nonlinear integro- differential 
equation, which accounts for both radiative and conductive heat transfer in the material. In the 
boundary condition at r = R, a mass transfer term due to evaporation will be included. Published 
solutions of Fourier’s equation concern linear or linearized forms of the boundary condition at 
r = R and thus are unreliable when the change of temperature of the radiating surface exceeds a 
certain small limit. 

DIFFERENTIAL EQUATIONS AND BOUNDARY CONDITIONS 

The problem outlined in the above section is governed by the following relations: 

1. Fourier’s equation: 


BI = 1 3 

3t a 2 r 2 3r 



for 0 < r < R, t > 0, 


where a 2 - p c/K , the thermal diffusivity (Reference 3). 
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2. The initial condition at time t = o : 


T(r,0) - T q - constant for 0 < r < R . 


( 2 ) 


3. The symmetry condition at the center r = 0 : 

BT(Q, t) 


dr 


0 for t > 0. 


4. The heat balance equation at the surface r =R: 

3T(R, t) _ 


-K- 


3r 


a Q e 0 T 4 (R,t) for t > 0, 


which equates heat radiated from the surface to heat conducted to the surface. 


(3) 

(4) 


According to the theory of the parabolic differential equation, the solution is determined by the 
six quantities T 0 , T(R, t) , R, t, o- 0 e 0 /K, and K /pc , which depend on the three fundamental units 
length, time, and temperature. The 77 -theorem of dimensional analysis predicts that the solution 
sought depends on 6-3 = 3 nondimensional power-products of the six quantities mentioned above. 

By use of the parameter f , defined by 


£ = 


Kpc 

w ( a o e o) 2 V 


(5) 


these expressions may be selected as follows: 


T(R, t) 


p , and D - 


H 

aR 


R T 0 3 


T o e o 


( 6 ) 


Characteristic features of the solution are shown by £ the parameter governing the time-scale: 

(a) £ depends on the product of K, p , and c, and (b) £ is very strongly affected by changes of the 
initial value T 0 of the absolute temperature. 


SOLUTION 

The second author will discuss in a later publication the existence and the uniqueness of 
solutions of the parabolic problem as defined by Equations 1-4. It is well known (see Reference 1) 
that substitution of the transformation 


T(r. 


t) = 


u(r, t) 

r 


+ T 


o 


into the differential equation (Equation 1) yields Fourier's equation 

3u(r, t) _ _1_ B 2 u(r, t) 

Bt a 2 Br 


(V) 

( 8 ) 
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for the planar case. The conditions in Equations 2, 3, and 4 are transformed as follows by use 
of Equation 7: 


u(r, 0) - 0 for 0 < r < R , 

Uo (t! 7 - 77 ) = 0 for t > 0 


( 9 ) 

(10) 


and 


for t > 0. 


f 1 3u(R, t) 
K \ R 3r 


n(R. t) \ 
R 2 / 


(t. * 


(u) 


Application of V Hospital T s rule reveals that the left-hand side in Equation 10 tends towards the 
same limit as 0.5 B 2 u/Br 2 as r - 0. This is equal to (a 2 /2)3u(0, t)/Bt, according to Equation 8. 
Therefore, the condition in Equation 10 may be replaced by 


u(0, t) = 0 for t > 0 , (12) 

By use of the complementary error integral (Reference 1), 


F(X) 



2 d/3 for 0 < X < oo and F(0) 


1 , 


(13) 


a function 


u ( r » t) 



V 


drj 



a(R - r ) \ dr i(v) 
2ft - 7) ) dr) 


dr) for t > 0, 


0 < r < R 


(14) 


is defined which depends on the free functions r 0 ( t) and Tj(t). These functions are assumed to be 
continuous and continuously differentiable, and may be interpreted as heat sources (or sinks) at 
r = 0 and at r = R, respectively. The right-hand side of Equation 14 satisfies Equations 8 and 9, 
since F(r,t)is a particular solution of Equation 8. The functions r 0 (t) and r x ( t) are determined 
by substituting Equation 14 into Equations 11 and 12 and by specifying that r 0 (0) = 0) = 0. 

Nondimensional variables are defined by 

t = fr* T o( t ) = T o and T i( t ) = T o R [“1 + 0^(0-)] • (15) 

By use of Equations 11, 12, 14, and 15, the following system of nonlinear integro- differential 
equations for cj 0 (<t) and aq(cr) is obtained: 

f" J 1 \ 

+ Jo F \2Dy5^F/ dr dr - 0 for o- > 0 , (16) 
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/ 


and 


r 

<H( T ) eXP ( 4 D 2 (cr -t)) 

d“i( T ) l 

Jo 

dr icr-T 

dr ya ~r 



F 


(■ 


2DVa 


i =) 


d"o ( T ) 

dr 


-j 4 

dr + ojjCcr) 


( 17 ) 


- D /n 


OJjCcr) 


-1 + 



for cr > 0 , 


where r is a dummy variable for o-. - Equations 14-17 confirm the above result of dimensional 
analysis. 

Differentiation of Equation 16 yields oj 0 ( 1 > (0) = 0 for i = 0, 1, 2, * • • and 0 < D < <»; i.e., the 
temperature changes initially only within an outer shell of the sphere under discussion. 

For D >0 ? the asymptotic expansion 


2 / 4 D\ 

^ 1 - + ■■■ 

is obtained. For D = 0 (i.e., R = °°),the right-hand side yields a series expansion 

Cl>j (cr) = 1 - (2/77 ■) y<T + ( 4/77)0- + ••• 


(18) 


(18a) 


in powers of |& , whose circle of convergence has been determined approximately. The range of 
applicability of Equation 18 depends on d; in case of D = o , this range is given approximately by 
0 < a < 0.05, according to numerical experience. 

Equation (17) is converted by use of the well-known equivalence of 

f(x) = 1 0i dZ and lo g(z)dz (18b) 


The equation resulting from this process and Equation 16 have been solved by use of the initial 
conditions oj 0 (0) = 0 and ^(0) = 1. The functions co 0 (a) and cj x (cr) were obtained at grid points 
cr > 0, i = 1, 2, * * * . For o < cr « l, the expansion in Equation 18 was used. The singularity at the 
upper limits r - a of the integrals was removed by use of an approximately valid closed-form 
evaluation of the integrals. 


4 



THE SEMI-INFINITE BODY (R^cd) 


Because of Equation 6, the limiting case R - °o is equivalent to D - 0 . By differentiating 
Equation 16 with respect to cr, it is seen that, for D - 0, d« 0 (o*)/dcr tends to zero uniformly in a- 
within every finite interval on the positive cr-axis. By application of the equivalence (18b), Equa- 
tion 17 becomes a nonlinear integral equation for ^(cr) : 

1 [ a °>i ( r ) 

i - ^iO) = ^rjo ~^== dr * for > °* wlth ^i^ 0 ) _ 1 • (19) 

Wolf and Mann (Reference 4) have shown both the existence and the uniqueness of solutions of the 
class of integral equations which includes Equation 19. According to this paper, the solution of 
Equation 19 satisfies the relations 


0 < ^(cr) < 1, (cr ) < 0, and ^(cr) - 0 


( 20 ) 


The solution of Equation 19 has been obtained numerically using the expansion in Equation 18a. 
A function co A (a) is defined as a solution of 


-if 

77 Jo 


“a 4 ( r ) 


dr - 0 for a > 0 


( 21 ) 


This equation is solved by 


u A (cr) - cr ” 1/8 for cr > 0 . (22) 

Because of Equations 20 and 21, o rl/8 may be interpreted as an approximately valid asymptotic 
solution of Equation 19 as cr - <». if co 1 (cr) in the left-hand side of Equation 19 as equated to a-’ 1/8 
is of the order of lCf <n+1) at cr = cr 0 > o, the function oj a (ct) = o-’ 1/8 satisfies Equation 19 with at least 
n decimals for a > cr 0 . 

NUMERICAL RESULTS 

The temperature T(r, t) can be expressed in terms of the constants T 0 , and D and by 
use of the numerical solutions for cj 0 (cr) and 60 1 (°0* According to Equations 5 and 6, £ and 
D depend on the initial temperature T 0 and geometrical and material properties which define 
the problem. In Figure 1, numerical results are presented for the time- dependence of tempera- 
tures at the surface and at certain stations inside spheres, with the radii 

R - 1 cm, R = 15.3 cm, and R ^ (23) 
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2000 


1000 



t, TIME (sec) 


Figure 1 — The curves with a vertical slope at t = 0 refer to surface temperatures; the other curves present tempera- 
tures at the center of the spheres or at the same distance from the surface of the semi -infinite body. The t-axis 
scale is linear over the range 0 < t < 1 and logarithmic for t >1 (sec). 


the initial temperature 


T 0 = 2000 °K 


(24) 


and the material properties 


K = 10 -4 kcal/m °K sec , 

p = 2000 kg/m 3 , 


c - 0.3 kcal Ag °K , 


e 0 = 0.8 , 


and 


cr 0 e 0 = 1.1 x 10 -11 kcal/(m 2 sec °K +4 ) 


(25) 
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The data of Equation 25 correspond approximately to the properties of quartz, which at T 0 = 2000 °K 
is still almost opaque to thermal radiation. The sphere with R = 15.3 cm weighs 30 kg. 

Figure 1 presents the surface temperatures for spheres with the three radii as given in Equation 
23, the temperatures at the centers of the two spheres with finite radii, the temperatures at the 
distances 1 cm and 15.3 cm underneath the surface of the semi-infinite body (R - a>),and the 
result of evaluating under the given conditions the approximately valid asymptotic solution 
(Equation 22) for the semi- infinite body. 

The slope BT(R, t )/Bt of the surface temperature tends to -a> as time t tends to zero. This 
singularity is a consequence of prescribing a constant initial temperature T 0 and a finite heat 
transfer rate as t - 0 (see Equations 2 and 4, respectively). In the solution derived, the singularity 
is caused by the complementary error integral (Equation 13) which exhibits different limits as 
r R for t = 0, and as t - 0 for r = R. . 

The temperature at a nonzero distance underneath the surface remains practically constant 
at the initial level T = T 0 for a time interval which increases together with this distance. Because 
of the converging cross-sectional areas as r decreases for finite R, the temperatures at the inner 
stations change more rapidly in the spheres of finite radius than in the semi- infinite body. With a 
small error increasing with time, the three surface temperatures presented in Figure 1 coincide 
until the temperatures at the centers of the spheres with finite radius R begin to drop. Initially, 
therefore, the semi- infinite solution furnishes a useful approximation for the surface temperatures 
of spheres with finite radius. In the semi- logarithmic presentation in Figure 1, the rates of 
decrease of all the temperatures presented are nearly linear with time for considerable intervals 
of time or temperature. An inspection of the graph suggests that the approximately valid asymp- 
totic solution (Equation 22) yields a reasonable continuation of the computed surface temperature 
of the semi -infinite body. 

The graph shows that the bodies made of quartz remain hot inside for a surprisingly long 
period of time. To check the order of magnitude of these results, a comparison to available 
linear solutions has been carried out by replacing Equation 4 by the linearized boundary condition 


BT(R, t) 
Br 


T 


3 

mean 


T(R, t) for t > 0 


(26) 


where T mean represents a suitably estimated constant number between the initial temperature 
T 0 and the surface temperature at the end of the time interval being considered. A heat trans- 
fer coefficient h is defined by 


h = °-o e o Tmean/K = constant . (27) 

The temperature history at the surface and at the center of spheres with finite radius R is presented 
as a function of T 0 , hR, and h 2 t/a 2 in Reference 2 and Reference 3, where a is defined in connec- 
tion with Equation 1. Comparison with Figure 1 of data taken from either of these sources shows 
approximate coincidence of the time elapsing before the temperature begins to change at the centers 
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of the spheres with finite radii, and the average rate of decrease of the temperatures at these 
centers. This agreement obviously depends on atrial and error choice of T mean in Equation 
26, and thus does not permit us to use linearized solutions in the problem under discussion except 
for sufficiently small time intervals. 

The numerical results presented pertain to constant material properties approximating those 
of quartz; these data show that internal heat conduction and emission of radiation from the surface 
furnish an inefficient way of cooling objects in space. If the cooling process begins at temperatures 
higher than 2000°K, internal radiative heat transfer will initially cause a rapid temperature decrease 
in the entire body; at a temperature close to 2000°K, however, internal radiation becomes insignif- 
icant in quartz, and the cooling process outlined in this paper takes over. 

In this paper, the only heat transfer considered inside the sphere with radius R is that by 
conduction. Greenland and Lovering (Reference 5) have considered the complementary case when 
heat conduction is absent and the material is sufficiently transparent to thermal radiation as to 
render T( r , t ) a function of time t only. From their heat balance equation 

- P cf 7T R 3 = 477R 2 ct 0 £ 0 T 4 (t) , ( 2 8) 

Greenland and Lovering obtain by straightforward integration the elementary closed-form 
solution 


1 - T (°)" 3 ] • (29) 

Because of its independence of the radial coordinate, T(r, t ) as given by Equation 29 may be 
used to determine a temperature history prior to applying the solution derived in this paper, which 
rests on the assumption of a constant initial temperature T(r, 0) = T 0 = constant. 

Goddard Space Flight Center 

National Aeronautics and Space Administration 
Greenbelt, Maryland, January 10, 1967 
188-45-01-01-51 
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